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A new phase field model of microstructural evolution is presented that includes the effects of elastic strain 
energy. The model's thin interface behavior is investigated by mapping it onto a recent model developed by 
Echebarria et al (J. Exploiting this thin interface analysis the growth of solid state dendrites are simulated 
with diffuse interfaces and the phase field and mechanical equilibrium equations are solved in real space on 
an adaptive mesh. A morphological competition between surface energy anisotropy and elastic anisotropy is 
examined. Two dimensional simulations are reported that show that solid state dendritic structures undergo 
a transition from a surface dominated [10] growth direction to an elastically driven [11] growth direction by 
changes in the elastic anisotropy, the surface anisotropy and the supersaturation. Using the curvature and strain 
corrections to the equilibrium interfacial composition and linear stability theory for isotropic precipitates as 
calculated by Mullins and Sekerka, the dominant growth morphology is predicted. 

PACS numbers: 



I. INTRODUCTION 

Recent years have seen the increased use of phase 
field modeling of free-boundary problems related to 
microstructure growth. This method has the advan- 
tage in that it avoids explicit front tracking by making 
phase boundaries spatially diffuse through the use of 
one or more order parameters, </>, which vary continu- 
ously across interfaces. This method of simulation has 
been very successful in modeling dendritic solidification 
dH |26| . In particular, Phase Field theory has confirmed 
the predictions of the microscopic solvability theory of 
dendrite growth That is, the anisotropy in the 

surface energy (7), although small, controls the growth 
rate and morphology of the resultant precipitate [20]. 

More recently the phase field approach has been used 
to study dendrite growth directions in cases where mul- 
tiple sources of anisotropy can compete to control the 
resultant morphological structure. Haxhimali et al. l30ll 
showed that the dendrite growth direction in Al-Zn al- 
loys can change continuously from (100) to (110) with 
the addition of solute. The authors explained their ob- 
servations as an interplay between the four fold and six 
fold anisotropy parameters lfl2Tl . was conjectured to vary 
as a function of solute. Analogously, Provatas et al. 12711 
modeled the change from (100) directed dendirites to a 
seaweed structure as the surface energy anisotropy com- 
peted with the direction of the temperature gradient. 

The phase field method has also been used to model 
the effects of elastic and plastic effects on microstruc- 
tural evolution |H [H H [lHH HI S . In these stud- 
ies the elastic strain energy was added to the phase field 
energy functional and the elastic modulus tensor (Cijki) 
was assumed to be a function of the phase field <f>. Many 
of these studies have focused on the interaction of elas- 
tic fields on grow th and coarsening rates of solid state 
particles lfl6l |3 3fl - These studies have shown that elas- 



tic effects play an important role in controlling growth 
rates and morphology, consistent with previous analyti- 
cal predictions of Laraia et al. IU4I1 and Voorhees et al 
ll35ll . Energy minimization techniques also showed that 
the elastic free energy drives the interface to become un- 
stable for mismatched elastic coefficients under both ap- 
plied and self strain lfl8i [l9ll . Other studies have also 
incorporated the effect of mobile dislocations in phase 
field modeling l23h . They have examined such effects 
as the role of dislocations on the coarsening rate during 
spinoidal decomposition ||24ll . 

Another phenomenon of both scientific and practical 
interest, and that has has not received enough theoreti- 
cal study, concerns the morphology and growth rate of 
solid state dendrites when two sources of anisotropy can 
compete. Anisotorpy can result from the surface energy 
(7) or the more interesting result of anisotropy emerging 
in the elastic modulus tensor (Cijki). Experimentally it 
has been found that dendritic structures can appear in the 
solid state |[l5[ l34ll provided specific elastic and kinetic 
relationships between the parent and precipitate phases 
are met. As discussed by Husain et al. 134]], solid state 
dendrites are made possible if the precipitating phase and 
the parent phase have similar atomic lattice structures, 
there is a high rate of atomic transfer between the two 
phases and the transformation rate is controlled by diffu- 
sion in the presence of an anisotropy. More exotic struc- 
tures can also be found through competing mechanisms 
leading to interface instabilities as shown by Yoo j37il . 

This paper studies solid state dendritic growth in the 
presence of competitive interactions between surface en- 
ergy and elastic anisotropy. It quantifies a transition be- 
tween dominant growth morphologies controlled by each 
of these two anisotropic effects. Section [Hi presents the 
phase field model used in this study. Section [HI] ex- 
amines the effect of isotropic strains on the equilibrium 
composition. Section [IV] provides an approximation to 
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the anisotropic strain field at the interface of the precipi- 
tate. Section[V]illustrates the effect of the surface energy 
anisotropy, elastic energy anisotropy and supersaturation 
on the growth morphology of dendritic growth. Finally, 
the transition between a surface energy controlled mor- 
phology and an elastic energy controlled morphology is 
quantified in Section[VH 



Z\ and Zq) is a function dependent on the elastic mod- 
ulus tensor. Explicit expressions for these prefactors are 
given in the appendix (AppendixlAli for the cubic elastic 
modulus tensor in two-dimensions. 

It is straightforward to show using the common tan- 
gent technique that the equilibrium composition, denoted 
Cf 9 , is modified by elasticity according to, 



n. PHASE FIELD MODEL 

This section derives a phase field free energy for solid 
state transformations, which assumed a dilute alloy free 
energy for the bulk chemical thermodynamics. To this 
is added an additional contribution to account for elas- 
tic free energy. Corresponding phase field equations are 
then derived from this free energy. It is shown that the 
phase parameters can be related according according to 
the thin interface asymptotics of Ref. J2l to a good ap- 
proximation. The thin interface limit allows for signif- 
icant efficiency of simulations of the phase field equa- 
tions, particularly when combined with novel adaptive 
mesh refinement algorithms JUSl- 

The bulk chemical free energy is given by 

G(4>, C, e«) - f(T m ) - S(4>)(T - T m ) 

+ ?!L( c \ nc - c ) + E{ ( t>) c + f el { ( j>) (1) 
v 

C, (j> and eij are the solute concentration, order param- 
eter and total strain fields, respectively. The parameter 
R is the gas constant, v is the molar volume, T is the 
temperature and T m i the melting temperature of compo- 
nent A of a dilute binary alloy. The elastic free energy is 
defined by Hooke's law, given by 



fel{4>) = Cij k l{(j)){£i 
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where C h e? = T ^^f^ is the equilibrium coexistence 
line corresponding to the parent phase in the absence of 



(5) 



_ 1v [Z3 + Z1) 



is the correction to the 



elasticity, G e/e = ^ 

phase diagram due to a local change in the strain and 

c eq 

k = -rwq is the partition coefficient. 



A. A Phase Field Model For Elastically Influenced Phase 
Transformations 



By applying the usual dissipative dynamics for the 
order parameter, mass conservation for concentration 
and strain relaxation for the strains, the following ki- 
netic equations for the order parameter <fr, undercooling 
U = e \Zk (where k is the partition coefficient and u is 
defined in below) and the strain fields e, j are derived: 
Phase Mobility: 



rA(h) 2 ^j- = V- [W 2 A{n) 2 V<j)} 



-V.(|V0| 2 ^M(n)M§) 
V d(V(f>) J 
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-4>)-\{l + B)U{l-<l) 2 ) 



(6) 



Chemical Diffusion: 



where e*J4>) is the phase dependent eigenstrain, given 



by 



i+gQ) r *A 

2 t ij 



1-3(0) *B 
2 ij 



Here e' 
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e*Sij is the hydrostatic lattice eigenstrain with e* = 
aA ~" B and a,i is the lattice parameters in each phase. 
The phase dependent elastic tensor is interpolated across 
the interface by C m = ^f±C? jh} + ^C^. 
A convenient choice for the interpolation function g((j>) 
which maintains the bulk phases at (j> = ±1 is given as 
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The elastic energy can be re-cast as a function of <?(</>), 
leaving the elastic free energy in the form of 

U = z 3 {g{4>)f + Z 2 (g (0)) 2 + Zt <?(<£) + Z (4) 
Each pre-factor in the polynomial of g(cf>) (i.e., Z3, Z2, 



^ = V • {Dq(4>)C b VU) 



+(1 + (1 - k)U)C b 1 + 



OlZn 



.(i + (i-ft)tO 



fc + i 

, 1 — A; 

Strain relaxation: 

day = d 5G(c/>,C,eij) 
dx 3 ~ dxj Sef. 



C b (l-k)J dt 
dC b 



The explicit form of B in equation|6]is 

aZ\ \ olZi 



B = 



- 1 



dt 



= 



(1 - k)C e b q ) 2C 6 (1 - k) 
aZi 



-M4) 



2Ct\l-k) 
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where Cb is the equilibrium composition corrected for 
the strain by equation^ i.e., Cb = C^ q — | n-k) w ^ 
C^ q being the equilibrium composition in the absence 
of elasticity a = ^? (Z 3 has been neglected here, as 
discussed below). The quantities W, r and A are the usual 
phase field parameters setting the dimensions of space 
and time JIH^]. The function A(h) = 1 + e 4 cos(46>) 
controls the four-fold surface energy anisotropy, where 6 
is the angle of orientation between the interface normal 
n and a reference axis. 

The dimensionless undercooling (u = *\Zk ) i n 
equation [7] is modified from its form in Ref. |2] to in- 
corporate an elastic correction to the equilibrium com- 
position. This results in 

e u = ^ (k + 1 - (1 - k)g(<f>) - - 5 (0) 2 )) I 

where Z 2 is the prefactor to g(<p) 2 in the elastic free en- 
ergy. The function modulates the diffusion through 
the interface correcting for the diffuse nature of <f> and is 
given by 

* = C b (fc + l-(l-fc)0-§^(l-0 2 )) (ID 

ZL-b 

Two-sided diffusion is controlled by the function 
= q(^) ^- where q(4>) modulates the diffusion in 
the two phases, to simulate equal diffusion coefficients 
q(4>) = 1. The dimensional diffusion coefficient is de- 
noted D.The value of elastic coefficients used herein is a 
actually reported as aCijki ■ 

The phase field equations above are derived in the limit 
where Z3 — > 0, the condition for a small difference in 
the elastic coefficients in either phase. This condition 
also holds for all materials in which ^ << 1 holds. 
Generally speaking Z\ > Z 3 even in the most extreme 
disparities of the elastic coefficients, so this assumption 
is not unreasonable. 

Following the procedure for the derivation of the phase 
field model presented in J2l, the constants W a , t and A 
are inter-related by the asymptotic analysis of Ref. l2(ill . 
This analysis maps the phase field model above, without 
elasticity, onto the sharp interface limit governed by the 
Gibbs-Thomson condition d n t = C eq — ACd (n)n — 
ACf3(n)V, where d is the capillary length, /3fe(n) is 
interface kinetic coefficient and AC is the concentra- 
tion jump across the to phase interface. Echebarria et al. 
showed that for the case of vanishing kinetic coefficient 
(appropriate for the study of solid state denderites) the 
following relations must be obeyed: |p- = 0.8839A and 
D = Dt/W 2 = 0.6267A where D is the diffusion coef- 
ficient in dimensional units (m 2 /s). These relationships 
arise by an expansion of the phase (<j>) and composition 
(C) inside the interface, which is matched to solutions 
outside the interface such as to emulate the sharp inter- 



face boundary conditions. 

It expected that for the phase field model presented in 
this work the thin interface relationships between W , 
t and A will be, to lowest order, the same as those 
for the model of Echebarria et al. based on two ob- 
servations. The first is drawn from the work of Yeon 
et al. 0], who showed that, to first order, the inclu- 
sion of elasticity is decoupled from the curvature and 
kinetic effects in the Gibbs-Thomson condition, which 
reads C int = C eq - ACd (n)n - AC(3 k (n)V - G eie . 
In their asymptotic matching procedure the derivation of 
d and (3^ in) remained unaffected by the presence of the 
G e i e contribution. 

The second observations validating the use of the thin 
interface analysis of Ref. J2l for the phase field model 
presented here is that for our simulations it was foind 
that ^ << 1 in Equation and B < 1 in equation 
The first quantity was found numerically to be at least 
an order of magnitude smaller than the other terms in 
its equation. It was likewise also found that B « 1 
(which includes Z2 <C 1). Both B and are zero in 
the limit of zero strain, but are also small for small lattice 
eigenstrains. For larger strains these terms are still valid 
as long as the interfacial velocity (V) is small. In effect 
a small velocity eliminates any excess kinetic effects that 
arise at the interface. In this work only small eigenstrains 
(e* w 0.005) and small growth rates (V < 1/j.m/s) are 
considered. 

It should be noted that the small variable B accounts 
for the effect of a variable strain field on the develop- 
ment of the precipitate. It plays an analogous role to the 
temperature correction of A that is used when simulating 
solidification of alloys with non-linear coexistence phase 
boundaries 0]. 



III. GROWTH OF AN ISOTROPIC SECOND PHASE 
PRECIPITATE WITH COHERENT INTERFACES IN 
AN ISOTROPIC PARENT PHASE 

This section examines the interface compositions of an 
isolated second phase precipitated into a parent phase, 
which is grown under isotropic conditions in both the 
parent and precipitate phases. The surface energy is 
made isotropic by setting the surface energy anisotropy 
coefficient 64 = 0. The elasticity equations are formu- 
lated in terms of cubic tensor coefficients, ie Cn,Ci2 
and C44 as described in appendix [A] For isotropic lin- 
ear elastic coefficients the cubic elastic terms are related 

by, 

C44 = ^(Cn - C 12 ) (12) 

where for this simulation C\\ = 1011 and C12 = 729 
(Cij is dimensionalized by a ) and the coherent hydro- 
static eigenstrain is set to e* = 0.005. The convergence 
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FIG. 1 : A circular isotropic coherent precipitate is grown in an 
alloy with initial composition of C = 0.04 with a hydrostatic 
lattice eigenstain of t* = 0.005. A cross section of the compo- 
sition is shown with the corresponding equilibrium composition 
corrected for strain, Cb (-j + The compositions are 

interpolated to the center of the interface and are found to be in 
good agreement with the sharp interface boundary condition. 



constant of the phase field equations is A = 3. The 
equilibrium composition is C^ q = 0.1 with an initial 
alloy composition of C = 0.04 and the solute parti- 
tion coefficient is set to a value of k = 0.1. The diffu- 
sion coefficients are set equal in both phases thus remov- 
ing the phase dependency in the diffusion coefficient (ie. 
q((f>) = 1) in equation|7] 

The total domain size simulated is 640014 7 on a side 
with periodic boundary conditions. The phase field and 
diffusion equations are solved using an explicit time step- 
ping algorithm on an adaptive mesh with a grid spacing 
of dx — 0.391 at the lowest level of refinement and a 
time step of dt = 0.01. The displacement field is solved 
by direct Gauss-Seidel iteration, which is found to re- 
quire (O(N)) operations per time step on an adaptive 
mesh 

The precipitate particle is grown and a cross section of 
the composition field solved by Equation [7] is shown in 
Figure [TJ This figure indicates the corresponding value 
of the equilibrium interfacial composition as calculated 
by equation|5] Cb is plotted as a function of the phase by 
interpolating it through the interface to its corresponding 

precipitate side value by Cb + ^^kj. The com- 

position field is interpolated to the interface described by 
the point where = both from inside the precipitate 
bulk and from outside the precipitate in the parent phase. 
The points are denoted C'^ and C e g and are found to 
have excellent agreement with the interpolations to the 
center of the interface. 



IV. APPROXIMATION OF THE STRAIN FIELD 
AROUND A PRECIPITATE WITH CUBIC ELASTIC 
COEFFICIENTS 

In this section the elastic field is analyzed around a 
circular precipitate where the cubic elastic coefficients 
are equal in both the precipitate and matrix phases. The 
anisotropy is entered into the cubic elastic coefficients 
by introducing a deviation from the isotropic relation in 
equation[T2]as defined here by, 

= \{C U -C 12 )+P (13) 

where (3 is the deviation from the isotropic elastic coef- 
ficients. Cix and C12 remain unchanged. While the an- 
alytical solution to the isotropic strain field has been de- 
rived for elliptical inclusions under a hydrostatic eigen- 
strain |@], a solution to an anisotropic precipitate under 
the same conditions is mathematically cumbersome do 
deal with ll29ll and is instead solved here numerically. 



Effect of Beta on Strain around a static particle size 




Angle 



FIG. 2: (Colour Online) The strain fields e xx , e yy and their 
sum are plotted for two particle radii with j3 = 200 in units 
of the model. The effect of the deviation from isotropy (/3 = 
0, 100, 200) on the trace of the strain tensor (e xx + e yy ) is il- 
lustrated by the solid lines. 

A circular precipitate in a parent phase with identical 
elastic coefficients in both phases is considered. The 
dimensionless elastic coefficients are set to values of 
Cix = 1011 and C12 = 729 (Cy is dimensionalized 
by a, ie. aCV,), the hydrostatic eigenstrain is set to 
e* = 0.005 and the concentration field is made constant 
(The elasticity here is not influenced by compositional 
effects and therefore any concentration field will pro- 
duce similar results). The deviation from elastic isotropy 
is studied by two controls, the particle radius and the 
strength of the deviation from isotropic elasticity {(J). 

The behaviour of the strain trace due to changes in 
the strength of the elastic anisotropy ((3) is studied by 
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FIG. 3: (Colour Online) Measured values of U amp plotted vs 
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holding the precipitate radius constant while increasing 
the value of (3. The amplitude of the trace of the strain 
tensor (e xx + e yy ) is found to have a strong sensitiv- 
ity to variations in (3. This is illustrated in Figure [2] for 
(3 = 0, 100, 200, where the individual strain components 
and their sum at the precipitate/matrix interface are plot- 
ted as a function of the angle 6, with zero representing 
the (10) direction. Notice the four fold symmetry of the 
trace of the strain. Although the strain energy depends 
strongly on the value of (3, the precipitate radius is found 
to have little effect on the amplitude of any perturbation 
to the strain trace. 



V. CONDITIONS INFLUENCING THE 
MORPHOLOGY OF PRECIPITATES 



This section characterizes three distinct controlling in- 
fluences on the selection of a dominant morphology of 
precipitated dendrites. These are the anisotropies of the 
elastic tensor and surface energy and the supersaturation. 
The effect of each of these parameters are systematically 
tested by increasing the strength of each parameter while 
holding the other parameters constant. 

In the following phase field simulations, A = 3 and, 
as required by the sharp interface analysis 10], the di- 
mensionless diffusion coefficient is set to D = 0.6237A. 
The temperature is set such that the equilibrium compo- 
sition is C^ q = 0.1 and the partition coefficient used 
is k = 0.1. The elastic coefficients were converted 
to units of the model by the elastic modulating factor 
a = 6.005 • 10~ 92 j- and in these units the elastic coeffi- 
cients are set to a value of Cu = 1011 and Cn = 729, 



where C,-, 



xC, 



Precipitate structures are grown in a system with 
periodic boundaries, where the system size is set to 
6400VFX6400W (W being the interface width, deter- 
mined along with r from the asymptotic analysis used). 
The precipitates grew to sizes of at most 2000I4 7 and 
the solution to the displacement field drops off as 1 /R. 
This justifies the claim that purely isolated precipitates 
are studied while using the periodic boundaries. The dif- 
fusion coefficients and elastic coefficients have no phase 
dependence. The grid spacing is set to dx = 0.391T4 7 
and the explicit time step is set to dt = O.Olr. 



The 9 dependence of the strain trace can be approxi- 
mated to lowest order by a single Fourier mode defined 

by, 



U A mpcos{40) 



(14) 



To measure the functional form for the amplitude of the 
strain trace (JJamp), the anisotropic strength ((3) is var- 
ied and the amplitude of the strain trace (Uamp) is mea- 
sured for the corresponding waveform at the interface 
of the precipitate. These measured values are fitted to 
a functional form, given by the equation, 



MP 



1 (3e* 



2C n + /3 



(15) 



The functional form of equation Q3] shows good agree- 
ment for all values of (3 in the regime where 
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Cn 



< 2 



as shown in figure [3] This range of anisotropies is well 
within the limits of the relative anisotropic strengths that 
are studied here. 



A. Elastic Anisotropy (J3) 



The elastic anisotropy emerges from the elastic ten- 
sor (C u , C 12 and C 44 = |(C n - C 12 ) + (3). The 
anisotropy of the tensor is varied by holding Cu, C\<x 
constant and varying C44 through changes in (3. Figure 
@] (a-c) shows the effect of increasing the strength of the 
elastic anisotropy by increasing C44 while holding sur- 
face energy anisotropy and supersaturation (controlled 
via the average composition C ) constant. From top to 
bottom the values of (3 used are (3 = 25, 100, 400 re- 
spectively. As can be seen in this figure, a small elastic 
anisotropy causes the surface energy to dominate and the 
dendrite grows in the [10] direction ( Figure|4](a) ). When 
the elastic anisotropy is increased to sufficient strength ( 
Figure |4](c) ) the dendrite grows in the [11] direction. 
When the anisotropies effectively destructively interact 
the resultant structure leads to an almost isotropic growth 
morphology ( Figure |4](b) ). 



(a) C = 0.04,/3 = 25 and e 4 = 0.01 
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(b) Co = 0.04,/? = 100 and e 4 = 0.01 




(c) Co = 0.04,13 = 400 and e 4 = 0.01 



FIG. 4: (Colour Online) Morphologies of growth precipitates 
vs. variation in the elastic anisotropy in the modulus tensor 
through an alteration in (3. The composition field is plotted on 
the adaptive mesh, blue being low concentration and yellow is 
a higher concentration. The total system size is 6400 W , plotted 
is an area of 2QQQW on one side. . 



(a) C = 0.05./3 = 400 and e 4 = 0.01 
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(b) Co = Q.05,/3 = 400 and e 4 = 0.03 




(c) Co = 0.05,/3 = 400 and e 4 = 0.05 

FIG. 5: (Colour Online) Morphologies of growth precipitates 
for a variation in the surface energy anisotropy by modification 
of e 4 . The composition field is plotted on the adaptive mesh, 
blue being low concentration and yellow is a higher concen- 
tration. The total system size is , plotted is an area of 
on one side. 



B. Surface Energy Anisotropy (£4) 

The surface energy anisotropy is entered into the 
model using the simple form for 4-fold surface energy 
7 = 7 (1 — e 4 cos(46*)). The effect of surface energy 
anisotropy is examined by varying e 4 and holding both 
C (or equivalently the super saturation VL = 0.555) 
and the elastic anisotropy ((3 = 400) constant. Fig- 
ure [5] (a-c) shows the effect of increasing the strength 
of the surface anisotropy under these conditions. From 



top to bottom in Figure [5] the values of the e 4 used are 
e 4 = 0.01, 0.03, 0.05 respectively. Analogously with the 
effect shown in section lVAl increasing the strength of e 4 
causes the morphology to transform from a preferential 
growth along the [11] direction ( Figure [5] (a) ) to that of 
the [10] direction ( Figure 0(c) )with a transition region 
where the precipitate structure is isotropic ( Figure |5](b) 
)■ 
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C. Supersaturation(fi) 




(a) Co = 0.04,/? = 200 and e 4 = 0.03 




(b) C = 0.06,(3 = 200 and e 4 = 0.03 
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c) Co = 0.08,13 = 200 and e 4 = 0.03 



FIG. 6: (Colour Online) Morphologies of growth precipitates 
for a variation in the supersaturation with both anisotropic 
strengths. e 4 and (3, held constant. The composition field is 
plotted on the adaptive mesh, blue being low concentration 
and yellow is a higher concentration. The total system size is 
640014^, plotted is an area of 2000W on one side. 

The supersaturation is varied by changing the initial 
(i.e. average) alloy composition C a . Figure [6] (a-c) 
shows the effect of decreasing the supersaturation while 
holding the strength of the surface anisotropy and elastic 
anisotropy constant. From top to bottom the values of the 
initial alloy composition used are C = 0.04, 0.06, 0.08 
respectively. For large super saturations (Figure [6] (a)) 
the growth direction is dominated along directions pre- 



ferred by the surface energy, i.e. the [10] directions. As 
the supersaturation is decreased a transition from the [10] 
growth direction to the [11] direction is observed ( Figure 
0(b) and (c) ). 



VI. CHARACTERIZATION OF THE 
MORPHOLOGICAL TRANSITION 

In this section we discuss a technique by which the 
dominant precipitate growth direction can be predicted. 
First, the point at which this transition occurs is defined 
and measured by examining the envelope of the precip- 
itate tips in R — 9 space. For a given value of the pre- 
cipitate radius, analysis of the amplitude of the envelope 
makes it possible to determine the critical surface en- 
ergy anisotropy (e^), for a specified elastic anisotropy 
(J3), where a morphological transition from [10] to [11] 
growth directions occurs. The elastic strain energy is, 
however, proportional to the precipitate area (in 2D) 
and the surface energy contribution varies as the parti- 
cle perimeter. As such the critical e£ is also a function 
of the precipitate size. The particle size dependency of 
el vs. j3 is then found by balancing the Gibbs-Thomson 
corrections corresponding to surface energy vs. elastic 
anisotropy. This critical radius scale is found to be pro- 
portional to the Mullins-Sekerka instability radius (R c ). 
Finally, a condition relating e^p as a function of (3 and R c 
at the transition is proposed. 



A. Defining the Transition Point 

The transition point that characterizes the controlling 
mechanism of growth morphology is defined as the point 
at which all competing anisotropics exactly cancel. Un- 
der this condition an isolated precipitate will grow (ide- 
ally) as a circle (a sphere in three dimensions) until the 
interface becomes unstable by the Mullins-Sekerka in- 
stability. While in this case the interface will become 
unstable, the envelope around the particle will continue 
to grow as a spheroid. It is this envelope that allows the 
point of transition between anisotropically controlled di- 
rections to be characterized. 

The concept of the precipitate envelope is illustrated in 
Figure [7] where the interface for 4 precipitates with val- 
ues of (3 = 100, 200, 300, 400, e 4 = 0.03 and il = 0.666 

C° q — C 

in R — 6 space are plotted, where Vt = . b ,,„°n . The 
dashed line shows the envelope surrounding the inter- 
face. As the magnitude of (3 approaches the transition 
point the amplitude of the envelope decreases, approach- 
ing zero. 
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FIG. 7: (Colour Online) R - 8 space of 4 precipitates with 
different (3 values, 64 = 0.03 and Supersaturation = 0.666. 
The dashed lines illustrate the corresponding envelope of the 
precipitate. As the control values approach a critical point, the 
envelope amplitude goes to 0. 



B. Measurement of the Critical Surface Energy 
Anisotropy - ef 

The critical surface energy anisotropic coefficient (de- 
noted ) is defined as the value of £4 at a given super- 
saturation (ft) and elastic anisotropy {(3) which results in 
an envelope amplitude of zero, is interpolated from 
the amplitudes of the precipitate envelopes obtained by 
varying £4 for given values of the supersaturation fl and 
elastic anisotropy f3. 

The envelope amplitude is approximated by measur- 
ing the difference of the total growth distance from the 
center of the precipitate along the [10] direction to the 
growth distance along the [11] direction. The transition 
point is the interpolated value for £4 such that these am- 
plitudes approach zero. Seven different supersaturations 
are considered here, where £4 is varied between 0.005 
and 0.05 and (3 is varied from 100 to 400. The enve- 
lope amplitudes are measured at arbitrary times, chosen 
in each case, however, such that the precipitate has out- 
grown any initial transients. 

For each value of (3 at each alloy composition C , the 
envelope amplitude is plotted vs £4. The inset to Figure[8] 
illustrates this for an alloy with an average composition 
of C = 0.06, and the deviation from elastic isotropy 
is characterized for values of (3 = 100,200,300,400. 
For each value of (3, the data is fitted linearly and is 
interpolated to the transition line to extract the critical 
surface anisotropic value £4 ((3, f2). A linear relation- 
ship between £4" and (3 is found for each of the super- 
saturations studied. These linear relationships are plot- 
ted in Figure [8] for the alloy compositions of C D = 
0.04, 0.05, 0.06, 0.07. FigureEthus predicts that 

cf = A C J (16) 
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FIG. 8: The relationship between and j3 is extracted by fit- 
ting their plots. Shown here are the vs (3 for 4 supersatura- 
tions, Q = 0.666, Q = 0.555, Q = 0.444, Q = 0.333. Inset: 
Interpolation of is extrapolated by fitting measured values 
for U amp -> Ofor Q = 0.444. 



6.84-10^ 5 ,A .o.5 = 8.33-10- 5 ,^ .06 = 10.68-10 ~ 5 and 
A0.07 = 16.98 • 10~ 5 for the supersaturations studied. 



C. Critical Tip Radii at the Transition Point 

The competing anisotropic effects controlling mor- 
phology cancel when all correction terms that are de- 
pendent on the interface normal angle (9) in the in- 
terfacial equilibrium composition ( C\ nt = C^ 9 — 
AC ca pM ar it y (1 + 15e4/(0)) — AC e j astic (1 - e e h(6)) 
) exactly cancel. Here AC ca pUiarity and AC e i as tic are 
the isotropic corrections to interfacial equilibrium com- 
position and £4 and e e represent the strength of the sur- 
face energy and elastic energy anisotropics. Assuming a 
linear fourier expansion with a 4-fold symmetry in both 
f(9) and h(8) (ie. cos (48)) the terms in the interface 
solute correction can be grouped by order of the fourier 
expansion. 

Cf, = C Q 9 — (AC ca piU ar ity + AC e 

lastic ) 

- (e S urf - tele) COs(49) (17) 

where e sur f = 15f4 AC ca piiiarityK and £ eie = 
£ e AC 'elastic are the relative anisotropic strengths of the 
surface energy correction and the elastic energy correc- 
tion respectively. The factor of 15 in the capillarity 
term comes from the stiffness of the capillarity, d (n) = 
d a (l + 15e4Cos(4#)). The elastic anisotropy strength, e e , 
is linked to the strength of the elastic anisotropy through 
(3 and is derived below (see Equationl2"Tl). 

For an isotropic morphology to emerge the coefficient 
of 003(6*) in EquationfTTlis required to vanish, ie. 



where the fitting parameter Ac a has values of A0.04 = 



tsurf ^ele — 



(18) 
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In EquationfTTlthe capillary term (AC ca piUary) contains 
a curvature correction, while the elastic term (AC e i as ti C ) 
does not. When these terms balance each other a curva- 
ture k is selected, which will be associated with a critical 



radius Rg ans (ie. n c 
next. 



1/Rf rans ). This is determined 



The capillarity correction is &C cap iiiarity = (1 - 
eq which is used to calculate e,„ r f as 



k)d C b 



-surf 

tsurf = 15C£(1 - fc)d e 4 K 



(19) 



The total elastic correction for cubic coefficients is 
defined by Equation |5j as AC e / astlc = = 



1 Cn+C 12 r *2 



(1 



j ~Y-k — c ^ ~~ — £* — L )- ^ n tne arjsence of elastic 
anisotropy the strain trace (e xx + e yy ) is zero and there- 
fore 



— ^' clastic , , ^ 

4 1 — k 



(20) 



e e ; e is calculated by consideration of the anisotropy in 
the strain field by substituting Equation[l4]for e xx + e yy , 
giving 



D. Calculation of the Transition Tip Radius using Linear 
Stability Approximation 



In the previous section a selected precipitate radius is 
derived based on the interpolated values for the critical 
surface energy anisotropy. However, this method relates 
e 4 to fi only once the value of Ac a is measured. A 
method to approximate the selected radius by theoreti- 
cal consideration of the Mullins-Sekerka linear stability 
analysis on an isotropic particle is now shown. 

Mullins and Sekerka in 1963 JMIHl predicted a crit- 
ical particle size i?£ =2 = HR* (R* = k is the in- 
stability mode and SI is the supersaturation) after which 
the particle interface becomes unstable. A particle at the 
transition point in our study can be considered to behave 
similarly to an isotropic particle in a supercooled matrix 
and the tip's radius of curvature is assumed to be pro- 
portional to this critical radius. Here the supersaturation 
of the precipitate under elastic strain is modified by the 
elasticity according to O e ; = , where C& is the 



_ C b -Cg 
(l-fc)C' 

equilibrium interface composition modified due to elas- 
ticity by Equation [5j This supersaturation is used in the 
linear stability analysis result to predict a minimum crit- 
ical radius of, 



Rms — 22—^- 



(25) 



1 uu 



1 - e e cos(A6) 



(21) 



where e e = and U amp is defined by Equation [TBI 

This results in e e ; e becoming 



1 Cn + C\2 a * 2 



8(1 -k) C u +P 



■fit* 



(22) 



Substituting Equations [T9l and |22] into Equation [T8l the 
critical radius of curvature required to maintain isotropic 
conditions is given by, 



R 



120(1 - k) 2 C l d 



(Cii + 0) 



F c 



~*2 



(Cn+Cia) (3 



(23) 



A fit to a selected critical radii is attained by substituting 
the fitted equation for the critical surface anisotropy co- 
efficient (Equation[T6|i and by choosing a reference point 
of fi = 0. This results in a relationship for the magnitude 
of the fj = transition radius that depends on concentra- 
tion through Aq given as, 



Rt 



120(1 -fc) 2 (C n )CX 



{C n +C 12 )e 



*2 



(24) 



The next subsection discusses a separate method by 
which this radius is estimated without the need for mea- 
sured values of Ac a ■ 



A comparison of this instability radius with the fitted 
radius of equation [24] shows excellent linear agreement 
as shown in figure[9] Equation [24] is calculated such that 




15 15 2C 25 30 



FIG. 9: The fitted critical tip radius for isotropic growth, from 
equation \24\ vs the onset Mullins-Sekerka wavelength for seven 
different supersaturations. 



C\\ + (3 — > C\i in Equation [23] The (3 dependence of 



R 



is introduced into Equation [25] through its rela- 



tionship to e e by the term Cl \+p - The final result for the 
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critical selected radius as a function of f3 and fl e i is 



R 



MS 
trans 



5 fl e l Cn 



(26) 



E. Calculation of the Critical Transition Point 

With the scale of the critical tip radius determined by 
linear stability theory as a function of supersaturation, 
the required measurement of the relationship between 
and (3 can be eliminated in Equation[23] This is done by 
substituting the linear stability prediction of the critical 
radius (Equation l26li into the equation for the selected 
transition radius (Equation [23}. The resulting relation- 
ship is solved for the critical surface energy anisotropy 
pre-factor (ef). This relationship is 
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Cl2 e* 2 3(l + -L) 

150(l-fc) 2 (Cn+/3)C^ ei Pl + C n j 



(27) 

Equation[27]defines a morphological transition line as 
a function of supersaturation (O), elastic anisotropy (J3) 
and the anisotropy of the capillarity (64). The transition 
lines for fl ei = 0.606, 0.479, 0.353, 0.226 are plotted in 
figure [TO] Precipitates grown above the transition line 
will grow in the [10] directions while growth for con- 
ditions below the line will grow along [11] directions. 
Some morphologies are overplotted above and below the 
transition line for SI = 0.606 in figure[TT]to further illus- 
trate the utility of Equationl27l 





FIG. 11: The critical transition line from equation \27\ Above 
the line the precipitates prefer to grow in directions that min- 
imize the surface energy [10], and below the directions which 
minimize the elastic energy [1 1 ]. The curve is plotted vs beta 
for Co = 0.04(Ti e ; = 0.606). Overplotted are 4 simula- 
tions corresponding to different values of £4 and (3. The plot- 
ted line represents the value of 64 , the transition point where 
anisotropics cancel out. Colour Online. 



VII. SUMMARY 

We have introduced a phase-field model for the study 
of the morphological development of elastically stressed 
solid state precipitates. We considered particles that have 
coherent interfaces and are under elastic self stress by a 
lattice mismatch eigenstrain. We used a new finite differ- 
ence based adaptive mesh refinement algorithm to solve 
the phase-field and strain relaxation equations thereby al- 
lowing for very rapid solution times. By investigating the 
effects of supersaturation, elastic anisotropy in the elas- 
tic tensor and anisotropy in the capillarity we developed 
a scaling relationship to predict which anisotropy will be 
dominant in the morphological evolution of the precipi- 
tate. It is interesting note is the effect of supersaturation 
on the selected precipitate morphology and growth direc- 
tions. 

We would like to thank the National Science and Engi- 
neering Research Council of Canada (NSERC) for finan- 
cial support and SHARC-NET for supercomputing sup- 
port. 



FIG. 10: The critical transition points from eauation \27\ Above 
the line the precipitates prefer to grow in directions that min- 
imize the surface energy [10], and below the directions which 
minimize the elastic energy [1 1 ]. Transition lines for different 
Q, e i are plotted. Colour Online. 



APPENDIX A: CUBIC ELASTIC FREE ENERGY 
COEFFICIENTS 

In the generalized elastic portion of the phase field free 
energy, as described by equation [4] f e i = Z 3 (g 3 (</>)) 3 + 
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%2 {93 i.4>) ) +%i93 (<j>) + Zo > several unknown terms Z3, 
Z 2 , Z\ and Zo are introduced. These coeffiecients are 
dependent on the particular values of the elastic modulus 
tensor in either of the precipitate or matrix phases. Pre- 
sented here are the explicit forms for these functions for 
two sided cubic modula and a hydrostatic elastic eigen- 
strain of the form, 



The second order coefficient to the elastic energy in 
terms of the phase is calculated to be 

Z2 = g(3(Clii + Clia) - (C2 n + C2 12 ))e* 2 
+-AC2U + C2 12 - Cln - Cl 12 )(e xx + e yy )e*(A3) 



e* 
e* 



and in the presence of equal elastic coefficients this term 
becomes a constant. 



The zeroth order component Zo) has no dependence on 
the phase at all and is calculated to be 



-(Cl 11 + C2 11 )(e x 



+ -(Clu + C2 11 )(e w - 

1 e* 
-^{Cl\2 + C2 12 )(e xx - —)(e yv 



2 ; 



+{Cl M + C2 44 )e 2 xy (Al) 

This pre-factor has no dependence of phase (nor con- 
centration) and therefore it does not appear in either the 
phase mobility equation (equation[6]) or the chemical dif- 
fusion equation (Equation |7J since the growth kinetics 
are dependent on differences of energy. It does however 
appear in the static elasticity equation (Equation[8]l. 

The first order component is the most prominent term 
in the model equations and is calculated to be 



Zy = i(3(ci u 
1 



- 2 (Cln 
1 



Cl 12 ) + C2 n + C2 12 )e* 
Cl 12 ){e xz 



- i (Cl u -C2 11 )(e 



+ e 2 ) 

yy> 



:(Cli 2 — C2 12 ) 



xx^-yy 



-(CI44 - C2a)4 y (A2) 



The third order component is calculated to be 



l(Cl n -C2 n + Cl 12 -C2 12 )e 



*2 



(A4) 



and has no dependence on the dynamic strain field. In the 
presence of equal elastic coefficients this term vanishes 
completely. 
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